explore_bbmsy

Exploring trends in B/Bmsy.

library(here)

## here() starts at /Users/frazier/github/ohiprep_v2021

library(tidyverse)

## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──

## ✓ ggplot2 3.3.5     ✓ purrr   0.3.4
## ✓ tibble  3.1.5     ✓ dplyr   1.0.7
## ✓ tidyr   1.1.4     ✓ stringr 1.4.0
## ✓ readr   2.0.2     ✓ forcats 0.5.1

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(broom)

key <- read_csv(here("globalprep/fis/v2021/output/fis_bbmsy_gf.csv"))

## Rows: 609501 Columns: 8

## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (3): stock_id, bmsy_data_source, RAM_gapfilled
## dbl (5): rgn_id, taxon_key, year, bbmsy, mean_catch

## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

scores <- read_csv(here("globalprep/fis/v2021/output/fis_bbmsy.csv"))

## Rows: 146602 Columns: 4

## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): stock_id
## dbl (3): rgn_id, year, bbmsy

## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

catch <- read_csv(here("globalprep/fis/v2021/output/stock_catch.csv")) 

## Rows: 90333 Columns: 8

## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (4): TaxonName, CommonName, Resilience, stock_id
## dbl (4): year, fao_rgn, TaxonKey, tons

## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

bbmsy <- left_join(scores, key) %>%
  filter(bmsy_data_source=="RAM") %>%
  filter(is.na(RAM_gapfilled))

## Joining, by = c("rgn_id", "stock_id", "year", "bbmsy")

head(bbmsy)

## # A tibble: 6 × 8
##   rgn_id stock_id            year bbmsy taxon_key bmsy_data_source RAM_gapfilled
##    <dbl> <chr>              <dbl> <dbl>     <dbl> <chr>            <chr>        
## 1      1 Istiompax_indica-…  2001  1.49    600217 RAM              <NA>         
## 2      1 Istiompax_indica-…  2002  1.47    600217 RAM              <NA>         
## 3      1 Istiompax_indica-…  2003  1.47    600217 RAM              <NA>         
## 4      1 Istiompax_indica-…  2004  1.44    600217 RAM              <NA>         
## 5      1 Istiompax_indica-…  2005  1.34    600217 RAM              <NA>         
## 6      1 Istiompax_indica-…  2006  1.35    600217 RAM              <NA>         
## # … with 1 more variable: mean_catch <dbl>

bbmsy <- bbmsy %>%
  select(stock_id, year, bbmsy) %>%
  unique() 

bbmsy_trends <- bbmsy %>%
  group_by(stock_id) %>%
  do(mdl = tidy(lm(bbmsy ~ year, data = .))) %>%
  unnest(mdl)  #fitting a model per each goal and rgn_id

## Warning in summary.lm(x): essentially perfect fit: summary may be unreliable

# data.frame(glance(data_lm, mdl))
bbmsy_trends2 <- bbmsy_trends %>%
  filter(term == "year") %>%
  select(stock_id, term, estimate, p.value) %>%
  data.frame() %>%
  filter(estimate <= 2) %>%
  filter(estimate >= -2)
hist(bbmsy_trends2$estimate)

stocks <- unique(bbmsy$stock_id)

for ( stock in stocks){ # stock <- stocks[1]
    f <- ggplot(filter(bbmsy, stock_id==stock), aes(x=year, y = bbmsy)) +
      geom_point() +
      geom_line() +
      xlim(c(2000, 2020)) +
      labs(title=stock) +
      theme_minimal() 
    print(f)
}

## geom_path: Each group consists of only one observation. Do you need to adjust
## the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to adjust
## the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to adjust
## the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to adjust
## the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to adjust
## the group aesthetic?